# Markdown
# library(bookdown)
# General data analysis and transformation
library(readr)
library(readxl)
library(dplyr)
library(janitor)
# library(gtools)
library(stringr)
library(units)
library(tidyr)
library(rlang)
# API-related
# library(jsonlite)
# library(urltools)
# Mapping and GIS operations
library(sf)
# library(leaflet)
# library(htmlwidgets)
# library(geojsonsf)
# library(rmapshaper)
library(tmap)
library(tmaptools)
library(ceramic)
# workflow
library(here)
# plotting
library(ggplot2)
library(forcats)
# Get list of CES 3 variable names and plain text names
ces_variables <- read_csv(here('data_processed',
'ces_names.csv'))
ces_variables <- ces_variables %>%
mutate(name_revised = make_clean_names(name))
ces_measures_overall <- ces_variables %>%
filter(group == 'overall',
subgroup == 'overall',
type == 'score')
ces_measures_effects <- ces_variables %>%
filter(group == 'pollution burden',
subgroup == 'environmental effects',
type == 'score')
ces_measures_exposures <- ces_variables %>%
filter(group == 'pollution burden',
subgroup == 'exposures',
type == 'score')
ces_measures_sensitive_pops <- ces_variables %>%
filter(group == 'population characteristics',
subgroup == 'sensitive populations',
type == 'score')
ces_measures_socioeconomic <- ces_variables %>%
filter(group == 'population characteristics',
subgroup == 'socioeconomic factors',
type == 'score')
This document attempts to do some quantitative geospatial analysis to look at correlations between redline mapping and indicators of water quality, public health, and facility information. It computes the area weighed average percentile for different CalEnviroScreen 3.0 (CES) indicators by HOLC rating for a given city and statewide. CES considers two broad categories of indicators, (1) Pollution Burden and (2) Population Characteristics, which are both further broken down into two sub-categories each (factsheet, indicators page).
The Pollution Burden category includes the following indicators:
The Population Characteristics category includes the following indicators:
# Define coordinate systems to use for transformations
projected_crs <- 3310 # see: https://epsg.io/3310
# load data
# Redline Polygons
redline_polygons <- st_read(here('data_processed', 'redline_polygons.gpkg'), quiet = TRUE)
# st_crs(redline_polygons) # to check the reference system
# transform to projected CRS
redline_polygons <- redline_polygons %>% st_transform(projected_crs)
# add the area of each redline polygon
redline_polygons <- redline_polygons %>% mutate(holc_poly_area = st_area(.))
redline_polygons <- redline_polygons %>%
mutate(holc_city = factor(holc_city),
holc_grade = factor(holc_grade))
# Redline cities list
holc_cities_list <- c('Fresno', 'Los Angeles', 'Oakland',
'Sacramento', 'San Diego', 'San Francisco',
'San Jose', 'Stockton')
# NOTE: not all redline polygons have an HOLC ID listed - create a unique ID for each one with a missing ID (so that each unique HOLC polygon can be associated with a weighted average CES score)
# create new column for the new holc id to be assigned
redline_polygons <- redline_polygons %>%
mutate('holc_id_2' = holc_id) %>%
arrange(holc_city, holc_grade, holc_id)
# initialiaze some variables for the loop below
grade_city_vec <- c()
j <- 1
# loop through all HOLC polygons ...
for (i_row in 1:nrow(redline_polygons)) {
cur_id <- redline_polygons$holc_id[i_row]
grade <- redline_polygons$holc_grade[i_row]
city <- redline_polygons$holc_city[i_row]
# the current combination of grade and city
grade_city_cur <- paste0(city, '_', grade)
# if the holc id is NA, check to see if there are already other missing holc id's for that city/grade combination - if so, increment a counter by 1, otherwise set counter (back) to 1
if (is.na(cur_id)) {
if(grade_city_cur %in% grade_city_vec) {
grade_city_counter <- grade_city_counter + 1
} else {
grade_city_counter <- 1
}
# enter a new unique ID for the row with the missing ID
redline_polygons$holc_id_2[i_row] <- paste0('unknown_',
grade, '_',
grade_city_counter)
# add the grade/city combo to a running list to track number of missing holc id's for that city/grade combination
grade_city_vec[j] <- grade_city_cur
j <- j + 1
}
}
# CES 3 Polygons
ces3_poly <- st_read(here('data_processed', 'ces3_poly.gpkg'), quiet = TRUE)
# st_crs(ces3_poly) # to check the reference system
# modify the 'City' column to better keep track when joining with other datasets
ces3_poly <- ces3_poly %>% rename('ces_city' = 'Nearby_City')
# clean names
ces3_poly <- ces3_poly %>% clean_names()
# Transform both to projected crs
# st_crs(redline_polygons) # redline_polygons <- redline_polygons %>% st_transform(crs = projected_crs)
# st_crs(ces3_poly) # ces3_poly <- ces3_poly %>% st_transform(crs = projected_crs)
# create a map to check that the data loaded correctly
# tm_shape(ces3_poly) + tm_borders() + tm_shape(redline_polygons) + tm_borders(col = 'red')
basemap <- function(city) {
# get a basemap (See: https://github.com/mtennekes/tmap/issues/185)
Sys.setenv(MAPBOX_API_KEY = Sys.getenv('mapbox_api_key'))
map_background <- cc_location(ces3_poly %>%
filter(ces_city == city),
verbose = FALSE)
map_background
}
holc_map <- function(city) {
background <- basemap(city)
# create map
map_redline <- tm_shape(background) +
tm_rgb(alpha = 0.5) +
tm_shape(redline_polygons %>%
filter(holc_city == city), is.master = TRUE) +
tm_polygons('holc_grade',
palette = c('green', 'blue', 'yellow', 'red'),
title = 'HOLC Grade') +
tm_layout(legend.bg.color = 'white',
legend.bg.alpha = 0.7,
legend.frame = TRUE)
map_redline
}
# Fix self-intersecting polygons
if (sum(!st_is_valid(ces3_poly)) > 0) {
ces3_poly <- sf::st_buffer(ces3_poly, dist = 0)
}
if (sum(!st_is_valid(redline_polygons)) > 0) {
redline_polygons <- sf::st_buffer(redline_polygons, dist = 0)
}
# get the intersection of the CES polygons and redline polygons
ces_redline_intersection <- st_intersection(ces3_poly, redline_polygons)
# inspect
# glimpse(ces_redline_intersection)
# add the area of each clipped polygon to the data frame
ces_redline_intersection <- ces_redline_intersection %>%
mutate(clipped_area = st_area(.))
# glimpse(ces_redline_intersection)
# # write to output file
# st_write(obj = ces_redline_intersection,
# here('analysis-output', 'ces_redline_intersection.shp'))
ces3_poly_map <- function(ces_id, ces_measure_name, city) {
background <- basemap(city)
# plot the CES polygons
map_ces <- tm_shape(background) +
tm_rgb(alpha = 0.5) +
tm_shape(ces3_poly %>%
filter(ces_city %in% (ces_redline_intersection %>%
filter(holc_city == city) %>%
pull(ces_city)))
) +
# filter(ces_city == 'Sacramento' | ces_city == 'West Sacramento')) +
tm_polygons(col = ces_id,
title = ces_measure_name,
border.alpha = 1) +
tm_shape(ces_redline_intersection %>%
filter(holc_city == city),
is.master = TRUE) +
tm_borders(lwd = 0, alpha = 0) +
# tm_shape(redline_polygons %>% filter(holc_city == 'Sacramento')) +
# tm_borders(lwd = 2, col = 'black') +
tm_shape(redline_polygons %>%
filter(holc_city == city, holc_grade == 'D')) +
tm_borders(lwd = 2, col = 'red') +
tm_shape(redline_polygons %>%
filter(holc_city == city, holc_grade == 'C')) +
tm_borders(lwd = 2, col = 'gold2') +
tm_shape(redline_polygons %>%
filter(holc_city == city, holc_grade == 'B')) +
tm_borders(lwd = 2, col = 'blue') +
tm_shape(redline_polygons %>%
filter(holc_city == city, holc_grade == 'A')) +
tm_borders(lwd = 2, col = 'green') +
# tm_text('holc_grade', size = 0.5) +
tm_layout(legend.bg.color = 'white', legend.bg.alpha = 0.7, legend.frame = TRUE)
map_ces
}
ces_redline_overlap_map <- function(ces_id, ces_measure_name, city) {
background <- basemap(city)
map_overlap <- tm_shape(background) +
tm_rgb(alpha = 0.5) +
# tm_shape(ces3_poly %>% filter(ces_city == 'Sacramento')) +
tm_shape(ces3_poly %>%
filter(ces_city %in% (ces_redline_intersection %>%
filter(holc_city == city) %>%
pull(ces_city)))
) +
tm_borders(lwd = 1) +
tm_shape(ces_redline_intersection %>% filter(holc_city == city), is.master = TRUE) +
tm_polygons(col = ces_id, title = ces_measure_name, border.alpha = 0) +
# tm_shape(redline_polygons %>% filter(holc_city == 'Sacramento')) +
# tm_borders(lwd = 2, col = 'black') +
tm_shape(redline_polygons %>%
filter(holc_city == city, holc_grade == 'D')) +
tm_borders(lwd = 2, col = 'red') +
tm_shape(redline_polygons %>%
filter(holc_city == city, holc_grade == 'C')) +
tm_borders(lwd = 2, col = 'gold2') +
tm_shape(redline_polygons %>%
filter(holc_city == city, holc_grade == 'B')) +
tm_borders(lwd = 2, col = 'blue') +
tm_shape(redline_polygons %>%
filter(holc_city == city, holc_grade == 'A')) +
tm_borders(lwd = 2, col = 'green') +
#tm_text('holc_grade', size = 0.5) +
tm_layout(legend.bg.color = 'white', legend.bg.alpha = 0.7, legend.frame = TRUE)
map_overlap
}
# FACET MAP
ces_redline_facet_map <- function(ces_id, ces_measure_name, city) {
background <- basemap(city)
map_facets <- tm_shape(ces3_poly %>%
filter(ces_city %in% (ces_redline_intersection %>%
filter(holc_city == city) %>%
pull(ces_city)))
) +
# tm_shape(ces3_poly %>% filter(ces_city == 'Sacramento')) +
tm_borders(lwd = 1) +
tm_shape(background) +
tm_rgb(alpha = 0.5) +
tm_shape(ces_redline_intersection %>%
filter(holc_city == city), is.master = TRUE) +
# tm_shape(ces_redline_intersection %>% filter(holc_city == 'Sacramento')) +
tm_facets(by = 'holc_grade', free.coords = FALSE) +
tm_polygons(col = ces_id, title = ces_measure_name, border.alpha = 0) +
tm_shape(redline_polygons %>% filter(holc_city == city)) +
tm_facets(by = 'holc_grade') +
tm_borders(lwd = 2, 'black') +
# # D
# tm_shape(redline_polygons %>%
# filter(holc_city == 'Sacramento', holc_grade == 'D')) +
# tm_facets(by = 'holc_grade') +
# tm_borders(lwd = 2, col = 'red') +
# # C
# tm_shape(redline_polygons %>%
# filter(holc_city == 'Sacramento', holc_grade == 'C')) +
# tm_borders(lwd = 2, col = 'gold2') +
# # tm_facets(by = 'holc_grade') +
# # B
# tm_shape(redline_polygons %>%
# filter(holc_city == 'Sacramento', holc_grade == 'B')) +
# tm_borders(lwd = 2, col = 'blue') +
# # tm_facets(by = 'holc_grade') +
# # A
# tm_shape(redline_polygons %>%
# filter(holc_city == 'Sacramento', holc_grade == 'A')) +
# tm_borders(lwd = 2, col = 'green') +
# tm_facets(by = 'holc_grade') +
tm_legend(bg.color= 'white', frame = TRUE)
map_facets
}
ces_scores_plot <- function(ces_id, ces_measure_name) {
# group by city and holc rating, then compute the weighted average score for each
ces_redline_grouped <- ces_redline_intersection %>%
st_drop_geometry() %>%
mutate(area_x_score = clipped_area * (!!as.name(ces_id))) %>%
group_by(holc_city, holc_grade) %>%
summarize(total_area = sum(clipped_area),
area_x_score_total = sum(area_x_score)) %>%
mutate(weighted_score = area_x_score_total / total_area) %>%
drop_units()
ces_redline_grouped_statewide <- ces_redline_grouped %>%
group_by(holc_grade) %>%
summarize(total_area = sum(total_area),
area_x_score_total = sum(area_x_score_total)) %>%
mutate(weighted_score = area_x_score_total / total_area) %>%
mutate(holc_city = '**Statewide**') %>%
drop_units()
ces_redline_grouped_city_state <- bind_rows(ces_redline_grouped, ces_redline_grouped_statewide)
# plot the results by city and holc rating, and include statewide summary
g_city_state <- ggplot(ces_redline_grouped_city_state) +
aes(x = fct_relevel(fct_reorder(holc_city, weighted_score),
c('**Statewide**'),
after = 0),
y = weighted_score,
fill = fct_rev(holc_grade)) +
geom_bar(stat = 'identity', position = 'dodge') +
geom_vline(xintercept = 1.52, size = 0.5, color = 'grey50') +
scale_fill_manual(values = rev(c('green', 'blue', 'yellow', 'red'))) +
guides(fill = guide_legend(reverse = TRUE)) +
labs(fill = 'HOLC Rating',
title = ces_measure_name,
x = NULL,
y = "Area Weighted Average Score By City") +
coord_flip()
g_city_state
}
ces_summary_plots <- function(ces_id, ces_measure_name) {
# # TROUBLESHOOTING
# ces_measures_effects_scores <- ces_variables %>%
# filter(group == 'pollution burden',
# subgroup == 'environmental effects',
# type == 'score')
# i_map <- 4
# ces_measure_name <- ces_measures_effects_scores$name[i_map]
# ces_id <- ces_measures_effects_scores$name_revised[i_map]
var_name <- ces_id # paste0('weighted_score_', ces_id)
# compute weighted average score for each HOLC ID (first compute for ces_redline_intersection polygon
holc_poly_wt_score <- ces_redline_intersection %>%
st_drop_geometry() %>%
mutate(area_x_score = clipped_area * (!!as.name(ces_id))) %>% # compute weighted avg score for each ces_redline_intersection polygon
group_by(holc_city, holc_id_2, holc_grade) %>% # compute weighted avg score for each HOLC polygon (combining individual ces_redline_intersection polygons that are within each HOLC polygon)
summarize(total_area = sum(clipped_area),
area_x_score_total = sum(area_x_score),
component_poly_count = n()) %>%
mutate(!!var_name := area_x_score_total / total_area) %>%
ungroup() %>%
drop_units()
# for each HOLC polygon, compute departure of weighted average score from the citywide average for the city it's in
holc_poly_wt_score_departures <- holc_poly_wt_score %>%
group_by(holc_city) %>%
mutate(!!paste0(var_name, '_city_avg') := mean(!!as.name(var_name))) %>%
# summarize(city_avg = mean(!!as.name(var_name))) %>%
ungroup() %>%
mutate(!!paste0(var_name, '_departure') := !!as.name(var_name) - !!as.name(paste0(var_name, '_city_avg'))) %>%
select(-total_area, -area_x_score_total, -component_poly_count) %>%
{.}
# compute z scores for each holc polygon, by city
holc_z_scores <- holc_poly_wt_score_departures %>%
group_by(holc_city) %>%
mutate(!!paste0(var_name, '_sd_city') := sd(!!as.name(paste0(var_name, '_departure')))) %>%
mutate(!!paste0(var_name, '_z_score_city') := !!as.name(paste0(var_name, '_departure')) / !!as.name(paste0(var_name, '_sd_city'))) %>%
ungroup() %>%
{.}
# # check - summary of average and median z scores by city
# holc_z_scores %>%
# group_by(holc_city) %>%
# summarise(med = round(median(cleanups_z_score_city), digits = 3),
# avg = round(mean(cleanups_z_score_city), digits = 3))
# # DELETE
# # box plot of z scores
# delete_boxplot <- ggplot(holc_z_scores,
# aes(x = holc_city,
# y = !!as.name(paste0(var_name, '_z_score_city')))) +
# geom_boxplot()
# delete_boxplot
# # for info purposes, compute the statewide average of each holc grade's departure
# statewide_avg_departures_bygrade <- holc_poly_wt_score_departures %>%
# group_by(holc_grade) %>%
# summarize(sw_avg_dep = mean(!!as.name(paste0(var_name, '_departure'))),
# sw_median_dep = mean(!!as.name(paste0(var_name, '_departure'))))
# box plot - raw scores
raw_scores_boxplot <- ggplot(data = holc_poly_wt_score_departures,
aes(x = holc_grade,
y = !!as.name(paste0(var_name)))) +
geom_boxplot(aes(fill = holc_grade), notch = TRUE) +
# scale_color_manual(values = alpha(c('green', 'blue', 'yellow', 'red'), 0.4)) +
scale_fill_manual(values = alpha(c('green', 'blue', 'yellow', 'red'), 0.6)) +
geom_jitter(color='black', size=0.6, alpha=0.5, width = 0.2) +
theme(legend.position="none") +
labs(x = 'HOLC Grade', y = paste0(ces_measure_name, ' Raw Score')) +
geom_blank()
# raw_scores_boxplot
# box plot - departures
departures_boxplot <- ggplot(data = holc_poly_wt_score_departures,
aes(x = holc_grade,
y = !!as.name(paste0(var_name, '_departure')))) +
geom_boxplot(aes(fill = holc_grade), notch = TRUE) +
# scale_color_manual(values = alpha(c('green', 'blue', 'yellow', 'red'), 0.4)) +
scale_fill_manual(values = alpha(c('green', 'blue', 'yellow', 'red'), 0.6)) +
geom_jitter(color='black', size=0.6, alpha=0.5, width = 0.2) +
theme(legend.position="none") +
labs(x = 'HOLC Grade', y = paste0(ces_measure_name, ' Departure')) +
geom_blank()
# departures_boxplot
# # density plot
# departures_density <- ggplot(holc_poly_wt_score_departures,
# aes(x = !!as.name(paste0(var_name, '_departure')))) +
# geom_density()
# departures_density
# point plot - all cities - departures
departures_point <- ggplot(holc_poly_wt_score_departures) +
geom_jitter(aes(x = !!as.name(paste0(var_name, '_departure')),
y = holc_city, color = holc_grade),
height = 0.25) +
scale_color_manual(values = alpha(c('green', 'blue', 'orange', 'red'), 0.3)) +
labs(x = paste0(ces_measure_name, ' Departure'), y = 'City') +
scale_y_discrete(limits = rev(levels(holc_poly_wt_score_departures$holc_city))) +
geom_blank()
# departures_point
# point plot - all cities - raw scores
raw_scores_point <- ggplot(holc_poly_wt_score_departures) +
geom_jitter(aes(x = !!as.name(var_name),
y = holc_city, color = holc_grade),
height = 0.25) +
scale_color_manual(values = alpha(c('green', 'blue', 'orange', 'red'), 0.3)) +
labs(x = paste0(ces_measure_name, ' Raw Score'), y = 'City') +
scale_y_discrete(limits = rev(levels(holc_poly_wt_score_departures$holc_city))) +
geom_blank()
# raw_scores_point
# point plot - departures - by city (facet on cities)
departures_point_city <- ggplot(holc_poly_wt_score_departures %>%
#filter(holc_city == 'Sacramento') %>%
{.}) +
geom_point(aes(x = !!as.name(paste0(var_name, '_departure')),
y = holc_grade,
color = holc_grade)) +
scale_color_manual(values = alpha(c('green', 'blue', 'orange', 'red'), 0.3)) +
scale_y_discrete(limits = rev(levels(holc_poly_wt_score_departures$holc_grade))) +
facet_wrap(vars(holc_city)) +
labs(x = paste0(ces_measure_name, ' Departure'), y = 'HOLC Grade') +
geom_blank()
# departures_point_city
# point plot - raw scores - by city (facet on cities)
raw_scores_point_city <- ggplot(holc_poly_wt_score_departures %>%
#filter(holc_city == 'Sacramento') %>%
{.}) +
geom_point(aes(x = !!as.name(paste0(var_name)),
y = holc_grade,
color = holc_grade)) +
scale_color_manual(values = alpha(c('green', 'blue', 'orange', 'red'), 0.6)) +
scale_y_discrete(limits = rev(levels(holc_poly_wt_score_departures$holc_grade))) +
facet_wrap(vars(holc_city)) +
labs(x = paste0(ces_measure_name, ' Raw Score'), y = 'HOLC Grade') +
geom_blank()
# raw_scores_point_city
# compute average departure for each HOLC grade within each city
holc_poly_wt_score_departures_summary <- holc_poly_wt_score_departures %>%
group_by(holc_city, holc_grade) %>%
summarize(city_holc_dep_average = mean(!!as.name(paste0(var_name, '_departure'))),
city_holc_dep_median = median(!!as.name(paste0(var_name, '_departure'))))
holc_city_means_plot <- ggplot(holc_poly_wt_score_departures_summary) +
geom_jitter(aes(x = city_holc_dep_average,
y = holc_city,
color = holc_grade),
height = 0.15) +
scale_color_manual(values = alpha(c('green', 'blue', 'orange', 'red'), 0.6), name = 'HOLC Grade') +
labs(x = paste0(ces_measure_name, ' Average Departure'), y = 'City') +
theme(legend.position = 'bottom') +
geom_blank()
# holc_city_means_plot
holc_city_medians_plot <- ggplot(holc_poly_wt_score_departures_summary) +
geom_point(aes(x = city_holc_dep_median, y = holc_city, color = holc_grade))
# holc_city_medians_plot
# return the plots from the function
list(departures_boxplot, departures_point, raw_scores_point,
departures_point_city, raw_scores_point_city, holc_city_means_plot,
raw_scores_boxplot)
}
i_city <- 1
ces_measure_name <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name) %>%
{.}
ces_id <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name_revised) %>%
{.}
map_redline <- holc_map(city = holc_cities_list[i_city])
map_ces <- ces3_poly_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
ces_redline_facet_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
i_city <- 2
ces_measure_name <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name) %>%
{.}
ces_id <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name_revised) %>%
{.}
map_redline <- holc_map(city = holc_cities_list[i_city])
map_ces <- ces3_poly_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
ces_redline_facet_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
i_city <- 3
ces_measure_name <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name) %>%
{.}
ces_id <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name_revised) %>%
{.}
map_redline <- holc_map(city = holc_cities_list[i_city])
map_ces <- ces3_poly_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
ces_redline_facet_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
i_city <- 4
ces_measure_name <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name) %>%
{.}
ces_id <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name_revised) %>%
{.}
map_redline <- holc_map(city = holc_cities_list[i_city])
map_ces <- ces3_poly_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
ces_redline_facet_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
i_city <- 5
ces_measure_name <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name) %>%
{.}
ces_id <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name_revised) %>%
{.}
map_redline <- holc_map(city = holc_cities_list[i_city])
map_ces <- ces3_poly_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
ces_redline_facet_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
i_city <- 6
ces_measure_name <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name) %>%
{.}
ces_id <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name_revised) %>%
{.}
map_redline <- holc_map(city = holc_cities_list[i_city])
map_ces <- ces3_poly_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
ces_redline_facet_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
i_city <- 7
ces_measure_name <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name) %>%
{.}
ces_id <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name_revised) %>%
{.}
map_redline <- holc_map(city = holc_cities_list[i_city])
map_ces <- ces3_poly_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
ces_redline_facet_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
i_city <- 8
ces_measure_name <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name) %>%
{.}
ces_id <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name_revised) %>%
{.}
map_redline <- holc_map(city = holc_cities_list[i_city])
map_ces <- ces3_poly_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
ces_redline_facet_map(ces_id, ces_measure_name, city = holc_cities_list[i_city])
ces_measure_name <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name) %>%
{.}
ces_id <- ces_variables %>%
filter(id == 'CIscore') %>%
pull(name_revised) %>%
{.}
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
# box plots
summary_plots[[7]]
summary_plots[[1]]
# point plots by city
summary_plots[[3]]
summary_plots[[2]]
# facet plots
summary_plots[[5]]
summary_plots[[4]]
# average scores by grade
summary_plots[[6]]
i_map <- 1
ces_measure_name <- ces_measures_effects$name[i_map]
ces_id <- ces_measures_effects$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 2
ces_measure_name <- ces_measures_effects$name[i_map]
ces_id <- ces_measures_effects$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 3
ces_measure_name <- ces_measures_effects$name[i_map]
ces_id <- ces_measures_effects$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 4
ces_measure_name <- ces_measures_effects$name[i_map]
ces_id <- ces_measures_effects$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 5
ces_measure_name <- ces_measures_effects$name[i_map]
ces_id <- ces_measures_effects$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 1
ces_measure_name <- ces_measures_exposures$name[i_map]
ces_id <- ces_measures_exposures$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 2
ces_measure_name <- ces_measures_exposures$name[i_map]
ces_id <- ces_measures_exposures$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 3
ces_measure_name <- ces_measures_exposures$name[i_map]
ces_id <- ces_measures_exposures$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 4
ces_measure_name <- ces_measures_exposures$name[i_map]
ces_id <- ces_measures_exposures$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 5
ces_measure_name <- ces_measures_exposures$name[i_map]
ces_id <- ces_measures_exposures$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 6
ces_measure_name <- ces_measures_exposures$name[i_map]
ces_id <- ces_measures_exposures$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 7
ces_measure_name <- ces_measures_exposures$name[i_map]
ces_id <- ces_measures_exposures$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 1
ces_measure_name <- ces_measures_sensitive_pops$name[i_map]
ces_id <- ces_measures_sensitive_pops$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 2
ces_measure_name <- ces_measures_sensitive_pops$name[i_map]
ces_id <- ces_measures_sensitive_pops$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 3
ces_measure_name <- ces_measures_sensitive_pops$name[i_map]
ces_id <- ces_measures_sensitive_pops$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 1
ces_measure_name <- ces_measures_socioeconomic$name[i_map]
ces_id <- ces_measures_socioeconomic$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 2
ces_measure_name <- ces_measures_socioeconomic$name[i_map]
ces_id <- ces_measures_socioeconomic$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 3
ces_measure_name <- ces_measures_socioeconomic$name[i_map]
ces_id <- ces_measures_socioeconomic$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 4
ces_measure_name <- ces_measures_socioeconomic$name[i_map]
ces_id <- ces_measures_socioeconomic$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
i_map <- 5
ces_measure_name <- ces_measures_socioeconomic$name[i_map]
ces_id <- ces_measures_socioeconomic$name_revised[i_map]
# map_ces <- ces3_poly_map(ces_id, ces_measure_name)
# map_overlap <- ces_redline_overlap_map(ces_id, ces_measure_name)
#
# tmap_arrange(map_redline, map_ces, map_overlap, nrow = 1)
#
# ces_redline_facet_map(ces_id, ces_measure_name)
#
# ces_scores_plot(ces_id, ces_measure_name)
summary_plots <- ces_summary_plots(ces_id, ces_measure_name)
summary_plots[[1]]
summary_plots[[4]]
# ces_measures_combined <- bind_rows(ces_measures_overall,
# ces_measures_effects,
# ces_measures_exposures,
# ces_measures_sensitive_pops,
# ces_measures_socioeconomic)
#
# for (i_output in seq(nrow(ces_measures_combined))) {
# ces_measure_name <- ces_measures_combined$name[i_output]
# ces_id <- ces_measures_combined$name_revised[i_output]
# var_name <- ces_id # paste0('weighted_score_', ces_id)
#
# # group by city and holc ID, then compute the weighted average score for each holc ID
# ces_redline_grouped_holc_id <- ces_redline_intersection %>%
# mutate(area_x_score = clipped_area * (!!as.name(ces_id))) %>%
# group_by(holc_city, holc_id_2, holc_grade) %>%
# summarize(total_area = sum(clipped_area),
# area_x_score_total = sum(area_x_score),
# component_poly_count = n()) %>%
# mutate(!!var_name := area_x_score_total / total_area) %>%
# drop_units()
#
# ces_redline_grouped_holc_id <- ces_redline_grouped_holc_id %>%
# select(holc_city, holc_id_2, holc_grade, !!var_name) %>%
# rename(city = holc_city)
#
# # add to output df
# if (i_output == 1) {
# df_output <- ces_redline_grouped_holc_id
# } else {
# df_output <- left_join(x = df_output,
# y = ces_redline_grouped_holc_id %>%
# st_drop_geometry(),
# by = c('city', 'holc_id_2', 'holc_grade'))
# }
# }
#
# df_output <- df_output %>% arrange(city, holc_grade, holc_id_2)
#
#
# # # write the output shapefile
# # st_write(obj = df_output,
# # here('analysis-output',
# # 'redline_ces_scores',
# # paste0('redline_ces_scores',
# # '_',
# # Sys.Date(),
# # '.shp'))
# # )
#
# # write the output to a geopackage file
# st_write(obj = df_output,
# here('analysis-output',
# 'redline_ces_scores',
# paste0('redline_ces_scores',
# '_',
# Sys.Date(),
# '.gpkg'))
# )